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We consider the properties of a one dimensional fluid of brownian inertial hard-core particles, 
whose microscopic dynamics is partially damped by a heat-bath. Direct interactions among the 
particles are represented as binary, instantaneous elastic collisions. Collisions with the heath bath 
are accounted for by a Fokker-Planck collision operator, whereas direct collisions among the particles 
are treated by a well known method of kinetic theory, the Revised Enskog Theory. By means of a time 
multiple time-scale method we derive the evolution equation for the average density. Remarkably, 
for large values of the friction parameter and/or of the mass of the particles we obtain the same 
equation as the one derived within the dynamic density functional theory (DDF). In addition, at 
moderate values of the friction constant, the present method allows to study the inertial effects not 
accounted for by DDF method. Finally, a numerical test of these corrections is provided. 

PACS numbers: 05.40.-a,61.20.Gy.05.10.Gg 

I. INTRODUCTION 

Over the last few years suspensions of interacting Brownian particles have been the subject of vivid theoretical 
interest due to new accurate experiments probing their properties at nanoscale, down to the effects of the correla- 
tion shells and layering structures in the density distribution. There is a great variety of systems and problems of 
fundamental and applied interest, including dense polymer solutions in good solvents^, the sedimentation of latex 
spheres, to produce materials with engineered optical gaps^, the design of micro-fluidic devices, to handle colloidal 
particles 3 ^*^ or the crowding effects in the cellular cytoplasm 6 .. Whenever the systems may be considered to be at 
thermodynamical equilibrium, the theoretical analysis of such structures may be efficiently done within the Density 
Functional formalism, with well tested approximations to include the effects of the repulsive and attractive interac- 
tions between the particles, although the inclusion of hydrodynamic (velocity dependent) interactions is still an open 
challenge. The theoretical study of the dynamical properties of colloidal particles suspended in a solution of lighter 
particles is a much harder problem, often studied through the Langevin approach to Brownian motioni^, with 
the lighter particles represented by a bath providing a damping force, with friction constant 7, and a thermalizing 
stochastic noise. 

Two levels of description, both based on the Fokker-Planck equation, can be employed to analyze Langevin model 
for Brownian motion. In the first, the so called Kramers equation 10 which governs the evolution of the joint probability 
distribution of position and velocity, one keeps track both of velocities and positions of the particles, whereas in the 
second, the Smoluchowski equationii, one considers only the evolution of the probability distribution of position. 
In fact, the velocity distribution relaxes in a time span of the order of the inverse of the damping constant toward 
its equilibrium form and afterward remains stationary, so that Kramers phase-space description becomes somehow 
redundant and one can restrict attention on the evolution of the spatial distribution, governed by Smoluchowski 
equation. However, the passage from the Kramers phase-space description to the Smoluchowski positional description 
requires the adiabatic elimination of the fast velocity variable. Even for the simplest case of ideal non-interacting 
particles, the correct procedure was understood only in the late seventies due to the work of Wilemskii 2 . and Titulaer—. 
In particular, Titulaer showed that a modified Smoluchowski equation can be derived from Kramers equation by 
means of a systematic 7 -1 expansion of the Chapman-Enskog type. He obtained the corrections to the standard 
Smoluchowski equation in terms of 7 for an arbitrary time independent external potential. More recently, Bocquet et 
ah 14 ' 15 gave a pedagogical discussion of such a derivation using the multiple time-scale metho d 16 : 17 . The corrections 
to Smoluchowski equation for large, but finite, values of 7 represent the effects of the underlying inertial dynamics, 
over the fully damped limit, in which at any time the velocity of a particle, averaged over the realization of the 
random noise, is proportional to the external potential force, (v(t)} ~ F(x(t)), with no inertial memory of the value 
of (v(t')) for t' < t. In the non-interacting case these corrections to the Smoluchowski equation produce a gradient of 
the external force, which determines a non uniform acceleration of the particle and renormalizes the effective diffusion 
constant. 

In this present paper, we are interested in the role played by the forces between the particles, in particular by 
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those having a short range character, such as the impulsive forces between hard spheres; which have been usually 
neglected in previous studies. We want to answer the question whether inertial effects matter, in the dynamics of 
a system of interacting colloidal particles, and which are the corrections to the Smoluchowski equation in that case. 
Interactions are expected to modify the motion of the particles by restricting their trajectories, by inducing different 
accelerations, and correlating their velocities and positions. In the case of an over-damped dynamics, i.e. when 
7 — ► oo, we presented^ a dynamic equation that governs the probability density of finding a particle in a given 
position. Starting from the Smoluchowski equation for the distribution function of the positions of N particles, we 
introduced a closure based on the assumption that the dynamical pair correlations could be approximated by those of 
a reference equilibrium system characterized by the same density profile as the non equilibrium system. The resulting 
self-consistent description for the average density was encoded in a deterministic Dynamical Density Functional (DDF) 
equation: 

9p(x,t) r ,\T-,,ti0F n i[p] 



dt 



= DV 



Vp(x,t)+p(x,t)V( : m »? +0V ext (x)) (1) 
dp{x,t) 



where (3J- n i[p\ is the non-ideal part of the free energy functional, f3V ex t(x) is the external potential, both in = 
(fcsT)" 1 unitsi^. The diffusion coefficient satisfies the Einstein relation D — where m is the mass of the 

colloidal particles, T the absolute temperature, fcs the Boltzmann constant and 7 the friction constant. Notice that 
the essence of any DDF approach is to set an approximate scheme in which the dynamic two-particle distribution 
function p2(x, x' , t), required to include the interactions effects in the time derivative of p(x,t), is taken as fully 
determined by the instantaneous value of that density distribution, as already done in earlier treatments like Enskog 2 - 
method and its revisions^. The exact time evolution of p(x,t) in interacting systems could only be obtained from 
the knowledge of the full previous history of the density distribution^, but from the practical point of view, the use 
of DDF approximations seems to be well supported by the comparison of their predictions with Brownian Dynamics 
Simulations^™! 

Can we extend such a description to the case of systems, where the dynamics is not over-damped? One would 
expect a richer dynamics as compared to the purely diffusive dynamics of eq.Q. Does the momentum of the particles 
play a role? The one particle phase-space distribution function, P[x, v,t], is the natural candidate to replace the 
density p{x,t) in this extended description. Of course, the Boltzmann equation for P[x,v,t], which predates all 
non-equilibrium kinetic equations, applies only to very dilute gases and does not incorporates the interaction with 
an heat-bath. We shall consider both these aspects and show that it is possible to derive eq.Q) as the leading term 
of a 7 _1 expansion, starting from the full inertial dynamics. The leading corrections are also obtained as the next 
terms in the expansion. In the present paper we investigate numerically and analytically the problem in the simplified 
version of a one-dimensional colloidal fluid driven by an heat bath at fixed temperature. Although it may appear that 
the one dimensional model employed is not of direct practical relevance, our motivation derives not only from the 
great simplification of the resulting algebra and computer codes, but also from recent experimental work for colloidal 
particles in very narrow channels^. 

An outline of this article is as follows: we open section II with a presentation of the microscopic model of inertial 
interacting particles subject to stochastic dynamics. We then introduce the evolution equation for the single particle 
phase-space distribution function obtained by combining the effect of dissipative collisions with the heat-bath, which 
gives rise to a Kramers-Fokker-Planck contribution, with the effect of inter-particle collisions, described by an Enskog 
collision term. At this stage, we separate the space dependence from the velocity dependence of the phase space 
distribution functions by using the eigenfunctions of the Fokker-Planck operator as basis functions. As a result of 
such a projection procedure we obtain an infinite non-linear system of coupled equations for the velocity moments 
of the phase distribution function. In section III by means of the method of multiple-scales we construct a uniform 
expansion in the inverse friction parameter and obtain the equation of evolution for the particle density. In section 
IV we explore the consequences of such an equation with a simple application and discuss its relation with the DDF 
equation. Finally in section V we draw the conclusions. 



II. ENSKOG-FOKKER-PLANCK EQUATION 



Let us consider a system of heavy particles suspended in a solution of lighter particles. Due to their smaller mass, the 
solvent particles perform rapid motions so that their influence on the heavy particles can be described by a stochastic 
force. As a result of such elimination of microscopic degrees of freedom one can represent the motion heavy particles 
by means of stochastic Langevin dynamics. Here we consider a system of N particles moving in one dimension, under 
the action of an external force f e (x) and interacting elastically with a pair potential energy U(x — x'). The equations 
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of motions are: 



i (2) 
-mTH + f B (xi) - £ ' r? i + &(*). ( 3 ) 



including the effects of the solvent with the linear friction coefficient 7, and the stochastic white noise with zero 
average and correlation 

(6 (*)& («)) = 2 1 mk B T8 lJ 5(t - s) , (4) 

T is the "heat-bath temperature" and (■) indicates the average over a statistical ensemble of realizations^*!^. The 
elimination in (J3J of the rapid bath variables £j (t) leads to the Fokker-Planck equation^, in terms of the probability 
distribution function, P(x,v,t) for the position and velocity variables, 
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d , f e (x) d 



dx m dv 



d T d 



—P(x,v,t)+v— + J -^ — P{x,v,t)=j —v + -— P(x,v,t)+C[x,v,t,P 2 }. (5) 



dv m dv 2 



The l.h.s is the Liouville operator for the ideal gas, under the external force f e {%), the first term in the r.h.s. 
represents the heat-bath as the standard Fokker-Planck collision operator, and the last term represents the effect of 
the interactions among the particles, as a generic collision operator, 

C[x,v,t,P 2 ] = ~ f dx' f dv' dU{x ~ X,) P 2 (x,v,x',v',t), (6) 
m ov J J ox 

This operator C satisfies the first equation of the Bogoliubov-Born-Green-Kirkwood-Yvon (BBGKY) hierarchy^, 
which connects the evolution of the n-particle distribution function, to the distribution function for (n + I) particles. 
For interacting particles, the evolution equation |J5J for the one-particle distribution function P(x,v,t) depends on 
the two-particle distribution P 2 (x, v, x', v', t), and some approximate closure is required to obtain a workable scheme. 

Whenever U(x-x') is a smooth function of the particle separation, like for the ultra-soft repulsive potentials used to 
model the steric repulsion between polymersi, or in the long range attractive interactions from dispersion or screened 
ionic forces, we may follow a mean-field approximation P 2 (x, v, x' , v ', t) ~ P(x, v, t)P(x' , v', t), which reduces (JSJ to a 
partial differential equation for the one-particle distributionSiLSi, and the effects of the particle interactions may be 
directly integrated over velocities, with the density distribution p(x,t) = J dvP(x,v,t); and included as a molecular 

field, f m (x, t) = — J dx' 9U<y Q~ x <■ p(x', t), to be added to f e (x) in the l.h.s. of as a self-consistent, p(x, t) dependent, 
force field. 

On the other hand, sharp repulsive contributions between the particles cannot be included as a molecular field, 
since they imply very strong correlations between the relative position [x — x') and the relative velocity (v — v') over 
the range of the repulsive force; so that P 2 (x, v, x' , v', t) goes sharply to zero when x — x' goes into the repulsive core. 
For hard-rod particles, of length a, there is an infinite force acting on an infinitesimal range around x — x' = ±ct, and 
the collision operator C[x, v, t, P 2 ] is exactly represented^ by the following operator 

KE[xi,Vi,t]= dx 2 / dv 2 Q{vi2s)(v 12 s) [S(x 12 - sa)b 12 - S(x 12 + sa)]P 2 [xi, v x , x 2 , v 2 , t], (7) 

s=±i^ 

where is the Heaviside function, vi 2 — (v± — v 2 ), x% 2 — (x\ — x 2 ) and bi 2 is the scattering operator defined for 
arbitrary function X(v\, v 2 ) by 

b 12 X( Vl ,v 2 ) = X(v[,v 2 ), (8) 

which for hard-rods swaps the velocities, b\ 2 X(vi,v 2 ) — X(v 2 ,vi), thus generating a correlation between relative po- 
sition and relative velocity. The representation J7J) formally integrates I© over the instant of collision, and substitutes 
the direct effect of the force by the change from the pre-collisional to the post-collisional velocities. 

A standard approximation of the collision term (JZj) is to assume that atoms are uncorrelated immediately prior 
to collision, which is the essence of Boltzmann's "Molecular chaos hypothesis" , but are correlated after they collide, 
because the collision itself generates correlations^. The revised Enskog theory (RET), developed by van Beijeren and 
Erns1p2i, truncates the infinite BBGKY hierarchy by factorizing 



P 2 (x! , vi, x 2 , v 2 , t) = g 2 [xi , x 2 ; p)P(x\ , v\ , t)P{x 2 ,v 2 ,t), 



(9) 
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The spatial pair distribution function, gi\x\, xi \ p], reflects the local positional correlations in the fluid. For particles 
with both short-range repulsions and long-range tails, the simplest approximation would be to split the generic collision 
operator K e in © into a molecular field representation of the soft interactions and an effective hard-rods description 
(|7|) of the core repulsion, following the usual treatment for equilibrium properties, which goes back to van der Waals, 
and is still the most used scheme within the Density Functional Formalism. 

In the case of one dimensional elastic hard-rods the RET provides the following expression for the collision integral: 

K E [xi,vi,t] = / dv 2 Q(v 12 s)(vi2s) x (10) 

S=±]V 

{g2[xi,x\ — sa; n]P[xi, v[,t]P[xi — scr, v' 2 ,t] — gi[xi,x\ + sa; n]P[xi,vi,t]P[xi + sa, V2,t]} 

Whereas in Enskog's formulation the pair correlation function at contact was assumed to be that of an equilibrium 
fluid evaluated at the local density at some point in between the colliding atoms, in the RET instead the contact 
value of <?2 is assumed: i) to be a non-local equilibrium functional of the local density, ii) to depend on time only 
through the density p[x, t) and iii) to have the same form as in a nonuniform equilibrium state whose density profile 
is p(x, t). Fortunately, in the case of a one dimensional hard-rod system the exact expression for the equilibrium pair 
correlation at contact is known given any arbitrary equilibrium density profile and reads^*: 

g2[a:±(7;p]= l-,(I±f) ■ 

The density dependence occurs entirely via the local packing fraction rj(x, t) = f^°/2 dx'p(x' ) t). 
At this stage it is convenient to introduce the following dimensionless variables: 

rEt2, v = —, x s * F = 7 ^ (12) 
a vt a vt 



F e (X).^#, F m (X,r)^^^ (13) 



P(X,V,t) = av T P(x,v,t), K(X,V,T)=a 2 K E {x,v,t) (14) 

where vt = \/ ksT/m. 

In situations where T >> 1 particles lose memory of their initial velocities after a time span which is of the order 
of the inverse of the friction coefficient 7 so that the velocity distribution soon becomes a Maxwellian. On the other 
hand, during the same interval the coordinates of the particles suffer a negligible change, as one can see comparing the 
product of the thermal velocity vt by 7 -1 with the typical molecular size a. In this limit the Smoluchowski description 
of a system of non interacting particles, which takes into account only the configurational degrees of freedom, turns 
out to be adequate. However, for intermediate values of T inertial effects may come into play. The question is how 
do we recover a description similar to that provided by the DDF approach starting from a phase description? On 
physical grounds one could directly neglect the inertial term in equation Q and consider only the evolution of the 
position distribution, as the DDF does, but such an approach does not give a clue on how the the inertial effects can 
modify the dynamics. 

Accordingly, Kramers' evolution equation for the phase space distribution function can be rewritten with the help 
of relations (|12I13|) and with the definition of effective field F(X, r) = F e (X) + F m (X, r) as: 

1 d p &,V,r) = LFpP{Xj T) _ l^p (Xi Vj T) 1 F{X: t)^P{X, V, t) + Ik(X, V, t) (15) 



having introduced the "Fokker-Planck" operator LfpP(X,V,t) = + V P(X,V,t), whose eigenfunctions 

H y {V)H u {V) = t^( — ^-Y~§v^ ex P(~ \V 2 ) have non positive integer eigenvalues v = 0,-1,-2,... Solutions of eq. 
(|15fl , where position and velocity dependence of the distribution function are separated, can be written as: 

00 

P(X,V,t) =Y,MX,t)H„(V). (16) 
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Moreover, by multiplying K(X, V, r) by ^H n (V) / H$(V) and integrating with respect to V, one represents the collision 
term as 



K{X,V,r) = J2C v (X,r)H v (V). 



(17) 



After substituting (JTBJ and ltT7|l into eq. 1(15)1 we find 



dr 
dX 



+ Tv<f> u (X,T) - C v (X,t) H„{V) + 



(18) 



- F(X)<p v (X )T ) H v+l {V)+v 



,(X,t) 



dX 



(*„,,, -l)JZ"„-i(V) = 



Finally, by equating the coefficients of the same basis functions, H v , we obtain an infinite hierarchy of equations which 
differs from standard Brinkman's expansion^ by the presence of collision terms. 



A. Physical interpretation of the expansion 



Before considering in detail the method of solution, we digress on the physical interpretation of our equations. By 
identifying </>q(X,t) with the dimcnsionless particle density, n — pa, 0i(A, r) with the momentum flow density, J v , 
4>2 = Ek — n/2 with the deviation from the thermalized value of the kinetic energy, Ek being the kinetic energy 
density, expressed in reduced units, we can rewrite the first three equations: 



dn(X,T) _ dJ v {X,r) 



Ot 



dX 



(19) 



dJv(X,T) 

dr 



-TJ V {X, t)+F(X, t)u(X,t) 



dE k (X,r) 
2 — +C 1 (X,r) 



(20) 



dE k (X,r) 
8t 



= -2T 



E k (X, T )--n(X, T ) 



dJ k (X, 
dX 



+ F(X,t)J v (X, t)+C 2 (X, t) 



(21) 



where the kinetic energy flow is defined as Jk = J dVV s P(X, V, t). 

Using the result derived in Appendix B we can express the coefficients C n (X, r) as divergences. First we introduce^ 
the kinetic pressure life = 2Ek and second identify the collisional contributions to the pressure and to the energy current 



We arrive at 



and 



dn c (x,r) 



dX 



-Ci(X,t) 



dJ c (X,r) 



dX 



-C 2 (X,r). 



dJ v (X,r) 



Or 



= -TJ V (X, t) + F(X, t)ji(X, t) - 



d{n k (X,T)+Tl c (X,T)} 
dX 



dE k (X, T ) 
dr 



= -2T 



E k (A, r) - -n{X, r)j +F(X, t)J v (X 7 t) - 



d[Jk(X,r) + J c (X, T )} 
dX 



(22) 



(23) 



(24) 



where the presence of the source term n(A, r)/2 maintains the fluid at constant temperature. Notice that properties 
(122(1 are consequences of the local conservation of momentum and energy during the collisional process. In one- 
dimensional clastic systems in addition to mass, impulse and energy all higher moments of the velocity distribution 
are conserved quantities under collisions, because K{X, V, r) is a divergence. 

If the hierarchy of moment equations is truncated, by supplementing the constitutive equations, one recovers the 
analogue of hydrodynamic equations with dissipation. We also remark that in a uniform bulk system the collisional 
contribution to the pressure coincides with the pressure excess over the ideal gas pressure, since II C (A, t) = n 2 /(l — n) 
, having used ea. (|B4|) and the contact value, g\ = 1/(1 — n), of the bulk pair correlation. 
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B. Exact solution for the free ideal gas. 



We illustrate the nature of the solutions by a simple example, namely the free-expansion of a system where the 
collisional terms and the molecular force field are dropped. Let us remark, that even in that simple case, the time 
evolution of the inhomogeneous ideal gas, is not well described by any simple truncation of the hierarchy, for instance 
setting $>z{X,t) = in order to obtain a closed system of equations for the first three weight functions. The exact 
eigenfunctions of Kramers' equation are known, and they can be expressed as infinite series of the form 



P M (X, V, t) = exp(-/xIY) exp 



6± 

r 



d_ 

dx 



1 



A_ d 
Tdx' 



(25) 



where A + and A_ are the raising and lowering operators on the FP velocity eigenfunctions, respectively,, A±H V (V) = 

Hu±i(V). The functions (^\x,t), which fully define PM(X,V,r) , may be any generic solutions of the diffusion 
equation 

d 2 



d_ 

dr 



(26) 



which produces the time dependence to be scaled as t\ = t/T. Therefore, for r> 1 there is a clear separation between 
the fast time dependence of the exponential decay exp{— fJ.Tr) and the slow dependence of the function (/^(X, r). 
The eigenfunction associated with /i = has the explicit form 



P (0) {X, V, t) = exp 



A+d_ 



^(X,r)= Ha{ V)^- H ^ V)d ^ 



(o) 



H 2 {V) d 



2A0) 



r 



dX 



2!r 2 dX 2 



(27) 



and represents a slowly decaying density inhomogeneity, 4>q(X, t{), with small (order 1/r, 1/r 2 , ...), slaved pertur- 
bations of momentum, energy, etc., whose shapes are given by the successive derivatives of the density distribution 
with respect to X. Similarly, the eigenfunction associated with fi = 1 has the explicit representation 



P (1) (X,V,t) = exp(-rY) 



'h(vu^ H 2 (V) H 3 (V) 

Hi(V) d 2 ^ 

r dx 2 




H (V)' 



dX 
dX 



(28) 



where the first line in the r.h.s. has the interpretation of a current inhomogeneity <f>Q (X,t/T), which slaves higher 
order (energy,...) perturbations with decreasing amplitudes (1/r,...), while the second line in the r.h.s. has the same 

structure of the P( \X,V,t) eigenfunction with amplitude 4>^ = T^ l dx4^\ an d both terms have the fast decay 
of the exponential pre- factor. The physical interpretation of such a combination is that an initially pure current 
fluctuation, described by Hi(V)<j>i(X, 0) would die very fast, as exp(— IY), but leaving behind a density fluctuation 

proportional to T~ 1 dx4 l< i ^ (X, 0), which would evolve with the slow time n. The particular combination in (|28|) is 
such that it completely cancels that remnant density fluctuations, i.e. it orthogonalizes (X, V, r) to P^ (X, V, r) , 
and leaves a purely fast decaying form. 

The structure of the higher order eigenvalues follows the same pattern, p( 2 \X, V,r) is an energy fluctuation, 
decaying as exp(— 2IY), but it has to contain diagonalizing terms proportional to r~ 1 P^ and to r~ 2 p(°\ to leave no 
slower remnant behind. With arbitrary choice of cp^ (X, 0), for v = 0, 1, 2, we may describe any initial distribution 
of the ideal gas, whose time evolution would be given by the superposition of the decaying modes. These "excited" 
fi > modes decay with a fast transient decay toward the only slowing decaying \x — mode, which contains t\ as the 
only relevant time scale. Such a separation between fast decaying exponential modes, and a slow diffusive mode should 
be much more generic than the particular realization in the free ideal gas. Indeed, it emanates from the structure of 
eq. (|15|l . where the heat-bath term is associated with the diagonal operator of the form TLpp — —Tv, which contains 
a null matrix element (y — 0), while the remaining elements are proportional to T. The non-diagonal contributions 
(given by the streaming terms for the ideal gas, and by collisions in general) are independent of T. In the limit r> 1, 
the generic structure of the eigenfunctions, reflects the properties of the eigenfunctions of the TLpp operator, with 
corrections of order 1/r, that is combinations of exponential decays exp(— vTt) and slow functions, evolving with 
T\ = r/r, or slower. Therefore, from an arbitrary initial condition, the system would have a fast transient decay 
toward a slow mode, made of a density distribution, accompanied of slaved current, energy, etc... fluctuations, with 
magnitude proportional to inverse powers of T. In the next section we work out the leading contributions of the 
collisions to that slow mode, taking into account their non-linear character generates slower than t\ times scales, as 
the slow reaction to a slowly changing external force F(X, t/T) would do. 
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III. MULTIPLE TIME-SCALE ANALYSIS 



How can we construct the equivalent eigenfunction representation of Kramers' equation for a system of interacting 
particles? The method is provided by the multiple time-scale analysis, that we shall discuss hereafter. The multiple 
time-scale method is designed to deal with non uniformities in systems with more than a time scale. It has been 
shown that a straightforward expansion of the Kramers equation in powers of the small parameter T _1 does not lead 
to a uniformly valid results^*. In order to obtain a uniformly valid expansion, instead, one makes use of the presence 
of two different time scales in the problem. The first scale, is fast and corresponds to the time interval necessary to 
the velocities of the particles to relax to configurations consistent with their thermal equilibrium value. The second 
time scale is much longer and corresponds to the time necessary to the positions of the particles to assume their 
equilibrium configurations. 

In the multiple time-scale analysis one determines the temporal evolution of the distribution function P(X, V, r) 
in the regime << 1, by means of a perturbative method. In order to construct the solution one replaces the 
single physical time scale, r, by a series of auxiliary time scales (tq, n, .., t„) which are related to the original variable 
by the relations r n = T~ n r. Also the original time-dependent function, P(X, V,t), is replaced by an auxiliary 
function, P a (A, V, r , t\, ..), which depends on the t„, which are treated as independent variables. Once the equations 
corresponding to the various orders have been determined, one returns to the original time variable and to the original 
distribution. 

One begins by replacing the time derivative with respect to r by a sum of partial derivatives: 



d _ d 1 d 1 d 



2 



(29) 



First, the auxiliary function, P a (X, V, tq, 7~i, ..) is expanded as a series of T 1 



P a (X, V, To, 71, T 2 , ..) - Yl fI P a S) ( X > V > T-0,71,7-2, ..). (30) 

Similarly, the collision operator is expanded as: 

C a (X,r) = J2^ C s , a (X,70,71,72,..). (31) 
s=0 

Next, each term P (l is projected over the functions H u : 

oo 

PW (X, V, r , n,..) =J2^(X,t ,t 1 ,T2,..)H„(V) (32) 

The term C SjQ , represents the contribution of order T~ s to C a (X, r): 

C s , a (X,r)= E g 2 (X,X + l)G% v ^(X,T) 1 p mu (X + 1,t) (33) 

l+rn=s fl,f 

- g 2 (X,X - l)G<Z^(X,T)ip mv (X - 1,t) 

One substitutes, now, the time derivative l|29ll and expressions (I30|) - (I32II into eq. p5|) and identifying terms of the 
same order in in the equations one obtains a hierarchy of relations between the amplitudes ip sv . The advantage 
of the method over the naive perturbation theory, is that secular divergences can be eliminated at each order of 
perturbation theory and thus uniform convergence is achieved. 

We show, now, how the method works. We substitute eas. (|29|) - (|3Uf) into eq. (|15|l and equating the coefficients of 
the same powers of T. To order T° one finds: 



L 



and concludes that only the amplitude ipQQ is non-zero. 
Next, we consider terms of order T _1 and write: 



fp [$^„.ff„]=0 (34) 



L 



FP 



^f°° Ho + DxiPqqHi — Cq.\Hi — C0.2H2 — Cq 3P3 (35) 

OT 
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having introduced, for notational convenience, Dx = (dx — F(X,t)). Following the method of reference^, the 
amplitudes with v = and s > are set equal to zero. Such a choice, although not unique is sufficient to eliminate 
secular terms, i.e. terms containing a dependence on the slow time To- By equating the coefficients multiplying the 
same H v we find that since: 

-d^ = ° (36) 
the amplitude f/'oo is not a function of To. Therefore, also the amplitude ipu, which is given by the relation 

ipu = -D x ipon + C ,i, (37) 

does not depend on t , being a functional of ipoo , both through the linear operator Dx and through the effective field 
Co i) whose explicit form is given in section IV. The remaining two amplitudes, instead, vanish because to order 
the self-consistent terms vanish, Co, 2 = and G),3 = 0: 

^12 - ^Co, 2 = 0, V13 = ^C ,3 = (38) 

In particular, the vanishing of Co, 2 is a consequence of the traceless form (for an elastic hard-rod system) of v (see 
appendix A). A similar property yields Co, 3 = 0. 
To order T~ 2 we obtain the equation: 



L 



FP 



1p2lH! + 1P22H2 + 1P23H3 

^plH x + + D x ^ n H 2 + dxi>xiH Q - C hl H x - d, 2 H 2 - C lt3 H 3 

OTq OTl 



(39) 



from which we obtain the conditions: 



and 



-dxipn, (40) 



dr 



-i>2i + C lA . (41) 



Notice that, since the l.h.s. of eq. (PHl does not depend on r , as discussed after eq. <|3*7j) . the r.h.s. must vanish. 
Utilizing eqs. (j37(l and 1)40(1 we write: 

^ = dx[Dxi>oo - C 0A ] (42) 
By carrying on the procedure to order T -3 we obtain: 



dtp, 



0(1 



dT2 



= -d x ip2i = -d x Cxx 



(43) 



where we have used eq. 11411 to eliminate -021 ■ 

For the sake of completeness we write the third order correction T~ 3 and find: 



on 



cV 3 



= -d 



x 



(d x F)(Dxipoo-C 0tl ) + C 2 ,i- 



dd 



0,1 



8t\ 



dxC\^ 



(44) 



The time derivative appearing in the r.h.s. can be expressed in terms of spatial derivatives of the order parameter 
ipoo using eq. (|42|l and therefore could be computed. 

As a check of the method we have re-obtained perturbatively the exact solution in the ideal gas case. Moreover, 
equation l|44|l reduces to the modified Smoluchowski diffusion equation in a potential obtained by Titulaer— , who 
showed that, in the case of independent particles in a parabolic potential, it coincides with the exact solution up to 
order T~ 3 . 
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In the following, we shall truncate the expansion to second order. Collecting together the various terms and 
employing eq. (|29[1 to eliminate the time variables tq, t±, t 2 and restore the original time variable r we obtain the 
evolution equation: 

^ = fd x [D x i/> 00 -C , 1 -~C 1 , 1 ] (45) 

Clearly, the evolution equation H45|l for the amplitude ipoo{X,T), representing the key result of the present paper, has 
to be supplemented with a prescription for Cq 1 and C\ \ which is given explicitly in the next section. These terms 
represent collisions and involve the density and current amplitude, ^>oo and ^n, respectively. However, the latter 
quantity can be expressed by means of equation l|45|l as a functional of ip 00 . In this manner expression i|37|) forms a 
closed equation for the density profile. 

It is worth to remark that, while in the original hierarchy eq. I|19|) the various amplitudes were independent fields, the 
solution obtained in this section, being, in fact, the generalization to interacting systems of the zeroth eigenfunction 
of Kramers equation, imposes a constraint on each of the v > components. We used this property as an internal 
check of the present extension to colliding particles. Employing the constraint provided by relations (|37JI and l|40|l into 
the first equations of the hierarchy (|19|l . we have verified that to order T~ 2 indeed the method provides a solution. 
In other words the solution even in the presence of collisions can be represented only by the eigenfunction associated 
with the less negative eigenvalue. 



IV. EVOLUTION EQUATION AND ITS DDF LIMIT 



Let us solve, in the case of interacting particles, the evolution equation l|45|l for the amplitude ipoo(X,T), which 
corresponds to the density fluctuation. The collisional contributions of orders and T -2 are, respectively: 



and 



Co,! - -lM*> T ) [^(X, X + l)lM* + hr)-g 2 {X,X- l)ijj m {X - 1, r)] } 



d,i= ^oo(X,t) g 2 (X,X + l)ip n (X + l,T)+g2(X,X-l)^ n (X-l,T) 

- ^xx{X, t) \g 2 (X, X + l)iP 0Q {X + l,T)+g 2 (X,X-l)^ 0Q (X-l,T) 
Jit L 



(46) 

(47) 
(48) 



where we have employed the matrix elements Gq = — 1, Gq ± = 2/y/n and Gi_o = — Co,i and relation ipn = 
—Dxipoo + Co,i to evaluate these expressions . Notice that the self-consistent interaction term Go,i depends only on 
the amplitude ^oo(^t) of the Hq(V) component. It describes the contribution to the effective restoring force when 
the the velocity distribution of the colliding particles is Maxwellian. The term G^i, instead, accounts for collisions 
between particles whose velocity deviates from the equilibrium thermal distribution. One may visualize, such a term 
by imagining the collision as occurring between a thermalized particle, i.e. a particle with zero average momentum, 
and a particle carrying momentum. Indeed, the Langevin dynamics leading to the standard DDF equation describes 
only collision between perfectly thermalized particles. This can be seen, by using eq. I|45|l and neglecting the term 
C± t i. One can recognize that the following equation 



d^oo(X,T) 
dr 



F(X,t)iP 00 (X,t) 



1 d ( d^oo(X,r) 
TdXl dX 

ii 00 (X, t) [g 2 (X, X + l)^oo (X + 1, r) - g 2 (X, X - l)Voo(* - 1, r)] } 



(49) 



represents the governing equation of the DDF method, expressed in dimensionless guise. 

It is also worth to comment the fact that the short-range and the long-range contributions to the dynamics, 
contained in C StV {X 1 r) and F(X, r) respectively, do not appear on equal footing. This state of affairs is encountered 
also when studying the equilibrium properties of liquids and was first recognized by van der Waals. In the present 
dynamical approach we see that the difference originates in the fact that in the hard core term, K(X, V, r), the velocity 
dependence of the distribution function P(X, V, r) does not factorize as in the molecular field term. The effect of 
hard-core collision depends not only on the amplitude of the Maxwellian component of the velocity distribution, but 
on the full velocity distribution. Therefore, as far as the system is not fully thermalized we observe a force which 
has not counterpart in equilibrium systems. However, as the system relaxes the term Gi^i tends to zero, because its 
amplitude depends on the current ipu(X, r) = — Dx4>oo(X, r) + Cq.\{X, t) which vanishes at equilibrium. 



10 



A. Linear analysis 



The following simple example may give an idea of the role of the corrections to the DDF equation. We compare the 
analytical results of our theory with those obtained by computer simulations of the model described by eq. J3J for 
an ensemble of N hard particles stochastically driven, in a periodic box of size Particle positions and velocities 
within two consecutive collisions are updated according to a second order discretization scheme for the dynamics 
cq. (J5J. Averages over 10 4 realizations of the noise were taken. 

We perform the analysis of the evolution of a small initial perturbation Ap(X,r)a = (iPqq(X,t) — poa) and show 
that while the DDF predicts that the relaxation depends only on the time scale t/T, hence is universal, the present 
theory leads to a violation of this scaling. 

In the limit of vanishing perturbations, each Fourier component evolves independently, and decays to zero expo- 
nentially^. 

The characteristic relaxation time can be ascertained by substituting in eq. 145|) the trial solution -0oo {X, r) = 
po<7 + A(t) cos(fcX) and keeping only linear terms. The resulting equation reads 



8A(t) 



a(K) 



A(r)+0(A 2 ), 



Ot r 

and has an exponential solution with a wave-length dependent decay time 

a(K) 1 K 2 



r r s(K) 

where S(K) is the hard-rod equilibrium structure factor 



(l-e(K)) 



S(K) 



sin(K) „. N2 sin 2 (f) 
1 + 2 Po a^ r ^ + 4( Po <j' 2 



K 



K 2 



(50) 



(51) 



(52) 



with p — po/(l — Po a )i being the pressure of the uniform ID fluid divided by kgT. In physical units the decay time 
associated with the wave- vector q = K/a reads Dq 2 (l — e)/S(qa), where for the self-diffusion coefficient, D, we use 
the result for isolated Brownian spheres, given by the Stokes-Einstein equation: 



D = 



(53) 



r] being the viscosity of the suspending fluid. The correction to the DDF result, e(K), appearing in ea. (|51|l . reads 



e(K) 



&Po<J 

i\7r 



(54) 



and depends on the ratio = 4^ between the Enskog collision frequency and the heat-bath characteristic 

frequency. For a uniform one dimensional hard rod system the collision frequency is 



v T k B T 

we = 2—=p = 2\ p 

Wir V mir 



(55) 



The prediction of the present theory for the variation of the relaxation time with respect to the wave- vector of the 
perturbation is shown in Fig.l and compared with the DDF result. 

Figure 2 illustrates the comparison between the theoretical predictions and numerical simulations in the case of an 
initial sinusoidal perturbation of period K — 2-k/IA and several values of the dimensionless friction T. Instead of the 
data collapse predicted by DDF we observe deviations in the short-time regime. Only for large value of T, i.e. in the 
over-damped limit, we recover universality. 

In particular, one observes a slower relaxation of density fluctuations. The larger the collision frequency, the slower 
the decay. In other words the theory predicts that, at fixed T, collisions render the relaxation process slower. What is 
the physical origin of this slowing down? One can think that a current of momentum can occur either via a particle 
displacement, i.e. a density change, or through collisional transfer. However, in the latter case the momentum can 
travel a distance a without paying any price to the frictional force —rwyv. Such a mechanism renders this relaxation 
"channel" slower. The two type of relaxation processes are sketched in Fig. 3. The first process dominates when the 
system is close to equilibrium when the velocity distribution is well described by a Maxwellian. In the second collision 
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process, which is more relevant far from equilibrium, the two distributions have the same temperature, but not the 
same momentum. 

Two facts are worth to mention: a) the correction has a kinetic origin as can be seen from the presence of the 
mass in the last member of eq. (|55f) . When to — > oo the correction vanishes, being the inertial effect negligible. On 
the contrary, in the case of over-damped dynamics, e(K) — > 0, the mass does not appear explicitly in the diffusion 
coefficient D and only geometrical factors such as S(K) play a role. Secondly, the correction increases as the particle 
size, a, increases. 



We have considered the non equilibrium colloidal dynamics of a system of hard rods of mass to driven by a uniform 
heat bath. The evolution depends on the non dimensional parameter T -1 , proportional to the time span occurring to 
the velocity distribution to reach its equilibrium value, and on the packing fraction. This evolution is described by a 
Kramers equation for the phase-space density P(X, V, r) supplemented by a collision term, treated within the Revised 
Enskog Theory. Since the momentum degrees of freedom equilibrate much faster than the positional degrees of freedom 
it is reasonable to look for a description which contains only the latter variables. By employing the multiple-time scale 
method we have performed the expansion of Kramers-Enskog equation and obtained a modified Smoluchowski- 
Enskog equation for the density field. We found that the collision term gives a non-local coupling between density, 
momentum and energy fluctuations. However, the density field slaves the remaining fields. To lowest order in 
the present method yields the same evolution equation for the density as the one obtained within the DDF approach. 
The present derivation does not require the existence of any equilibrium density functional, but is based on kinetic 
theory arguments. Therefore, it can be applied to generic non-equilibrium systems, where the RET closure of the 
evolution equation for the phase-space distribution is physically sound. However, containing as a key ingredient the 
same equilibrium pair correlation as the DDF, the matching between the two methods is not too surprising. 

As discussed by Archer and EvansAS if the thermal equilibration occurs mainly via the solvent the deviations from 
the DDF should be negligible. Nevertheless, for atomic fluids the harshly repulsive potential might concur appreciably 
to the relaxation process and lead to significant effects which are beyond the limits of the DDF approach. 

Besides reproducing known results the present derivation provides systematic corrections to the DDF equation 
accounting for the deviation of the velocity distribution from the Maxwcllian. Hence, it can describe situations very 
far from thermodynamic equilibrium or even situations where a steady, but non-equilibrium state exists. 

The present method quite naturally lends itself to the following future applications and extensions: a) hard core 
systems whose spatial dimensionality is larger than one, b) systems of particles experiencing inelastic collisions, such 
as granular gases, where free energy functional approaches are not applicable^! and the RET closure provides a valid 
alternative, c) systems having a non-uniform temperature profil e) 38 ! 39 where the standard isothermal DDF approach 
cannot be applied, d) inclusion of higher order corrections in the inverse friction expansion T _1 accounting for currents 
associated with higher moments of the velocity distribution. 
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APPENDIX A: COLLISION INTEGRALS 



We consider explicitly the first three coefficients C a featuring in the series expansion 1)170: 




(Al) 
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with /Uq = 1, /ii = V and \i 2 = (V 2 — l)/2. Using the definition of K given by ea. (|llfl one finds the expressions: 



g 2 (X,X-l] 




C a (X,r) = g 2 (X,X+l)i dV» a {V) 

rO pod 

duuP{X, V, t)P(X + 1,u + V,t)+ / duuP(X, u + V, t)P{X + 1, V, t) 

Jo 

dVfi a (V) 

duuP(X, u + V, t)P(X - 1, V, t) + / duuP(X, V, t)P(X -1,u + V, t)1 } 

After substituting expansion into eq. (A2) and integrating over velocities one obtains: 

C a (X, t) - g 2 (X, X + l)^2 G%M X , t)4>u(X + 1,t)- 92 (X,X-1)J2 G«^(X, t)^(X - 1, r) 



//.// 



where the matrix elements „ are given by 



Got 



dVfi a {V) / duuH /Ji (V)H„(u + V) + / duuH v (V)H^,(u + V) 



(A2) 



(A3) 



(A4) 



The integral Co(X, r) = is zero, as required by the conservation of the number of particles during a collision, and 
indeed all G° „ = vanish. The explicit form of the v for n, v = 0, 1, 2 and a = 1, 2 are given by the following 
matrices: 



-1 
1 





and 



G ,, ,, — 



1 _2_ 

'2 v^F 

1 I 

2 



--f= ± 

V7T 2 



APPENDIX B: A USEFUL IDENTITY 

We prove hereafter that the collision kernel K (X, V, r) and thus the expansion coefficients C n (X, r) can be expressed 
as divergences. To this purpose we employ the following identity^: 



f 1 d 

S(X,X + Y) - S(X -Y,X) = / dz — S(X — (1 — z)Y, X + zY) (Bl) 

Jo dz 

= 4v I dzS(X -(l-z)Y,X + zY) 

(jX j 

and identify 

/oo 
dV 2 (V 1 -V 2 ) (B2) 
-OO 

{e(Vi - v 2 )P(x, v u t)p(x + y,v 2i t) + e{v 2 - v x )p{x + y, v u t)p(x, v 2 , t )} 

and setting Y = 1 rewrite ea. (|ll|l with the help of ea. (|B2(l as: 
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K(X,Vx,t) =—gxj dzg 2 (X-(l-z)Y,X + zY) J dV 2 {V x -V 2 )x 
{0(Vi - V 2 )P(X - (1 - z)Y, Vl,t)P(X + zY, V 2 ,t) 
Q(V 2 - Vi)P(X + zY, V u t)P(X - (1 - z)Y, V 2 , r)} 



(B3) 
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FIGURE CAPTIONS 

Fig.l Exponential relaxation time corresponding to a sinusoidal density modulation of wave- vector K. The uniform 
background is p§o = 0.6. The dashed line represents the DDF results, while the continuous line represents the result 
of the present theory with T = 10 and the dots the non interacting system. 

Fig. 2 Initial stage of the temporal evolution of the amplitude decay of a sinusoidal modulation for different values 
of the dimcnsionlcss parameter T. The DDF theory would predict a collapse of the data when using the time scale 
r/r, whereas the present theory predicts a T dependence. The circles represent the numerical results relative to 
r = 5, square symbols refer to T = 7, diamonds to T = 10, upper triangles r = 20 and down triangles to T = 100. In 
the inset we report the same data using the original time scale r. The straight lines correspond to the values of the 
relaxation time predicted by the linearized solution. 

Fig. 3 Sketch of the two different collision processes contributing to the term Co,i (process I) and to the term C\.\ 
(process II). 
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